Magnetic phases and transitions of the two-species Bose-Hubbard model 
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A model of two-species bosons moving on the sites of a lattice is studied at nonzero temperature, 
focusing on magnetic order and superfluid-insulator transitions. Firstly, Landau theory is used 
to find the general structure of the phase diagram, and in particular to demonstrate the presence 
of first-order transitions and hysteresis in the vicinity of a multicritical point. Secondly, an ex- 
plicit thermodynamic phase diagram is calculated using an approach based on a field-theoretical 
description of the Bose-Hubbard model, which incorporates the crucial effects of particle-number 
fluctuations. The maximum transition temperature to a magnetically ordered Mott insulator is 
found to be limited by the presence of the superfluid phase. 



I. INTRODUCTION 

Experiments with cold atoms in optical lattices have 
provided a new window into the physics of Mott insula- 
tor and superfluid phases, and the phase transitions be- 
tween them.— An important recent development in 
this area has been the demonstration of superexchangc 
interactions in a spin mixture of bosonic atoms^ suggest- 
ing that the goal of simulating magnetic Hamiltonians is 
within reach. The primary obstacle to this breakthrough 
remains the low temperatures (or entropies^) that must 
be reached before magnetic ordering occurs. 

Perhaps the simplest possible model exhibiting the 
physics of magnetic ordering in itinerant systems 
is the Bose-Hubbard model^ with two species of 
bosons^ i 10 ' 11 ^ 2 In the limit of strong repulsion U and 
weak hopping t, it can be described by an effective model 
of localized pseudospins, using a perturbative expansion 
in t/U. The coupling between these moments is in- 
duced by virtual tunneling of particles between neigh- 
boring lattice sites^ giving a scale J ~ t 2 /U (see also 
Section ["II Bp . The critical temperature for magnetic or- 
dering is therefore strongly suppressed for t -C U, and, as 
in the analogous case of fermions^ increasing it requires 
hopping strengths t/U large enough that the effective 
local-moment description breaks down. 

To describe this regime, Altman et alJ^ introduced 
a method that builds upon mean-field theory, and in- 
cludes the particle-number fluctuations that are vital to 
magnetic ordering. They thereby determined the zero- 
temperature phase diagram, including Mott insulating 
phases with both ferromagnetic and antifcrromagnctic 
order, in agreement with perturbation theory^ and also a 
superfluid phase, which cannot be described by an expan- 
sion in t/U. In recent work, Soyler et alJ^ used quantum 
Monte Carlo simulations to find the zero-temperature 
phase diagram, confirming most of the features found 
by Altman et al., and also observing phases with simul- 
taneous superfluidity and lattice-symmetry breaking^ 

In the present work, we address the phase diagram at 
nonzero temperature, with the main focus on the various 
types of order that are possible, including superfluid and 
magnetic states. The main contributions are as follows: 
Firstly, Landau theory is used to understand the gen- 



eral form of the phase diagram; most significantly, it pre- 
dicts a broad region of first-order transitions and hystere- 
sis. We then introduce a framework for the studying the 
phase structure of the model that uses a field-theoretical 
approach based on a strong-coupling expansion of the 
Bose-Hubbard model. We use this to present an alterna- 
tive calculation of the zero-temperature phase diagram, 
which is mathematically equivalent to that of Altman 
et al.^ but makes use of a quite different formalism. 
Finally, we demonstrate the extension to nonzero tem- 
peratures, by applying it to the calculation of the phase 
diagram for temperature T > 0. 

The approximation method that we use is based on 
the standard mean-field theory for the Bose-Hubbard 
modeli^ This can be derived by using a Hubbard- 
Stratonovich transformation^ to write the partition 
function as an integral over an auxiliary field the 
mean-field theory is given by a saddle-point approxima- 
tion for this integral. In the Mott insulator, the integral 
is peaked at ip = 0, and the phase transition to the su- 
perfluid phase is signaled by a change to a nonzero ex- 
pectation value of tp and hence superfluid order. 

While this mean-field theory correctly predicts the 
phase structure of the spinless Bose-Hubbard model, it 
is unable to distinguish different spin-orderings within 
the Mott insulator. As noted by Altman et al.,— this 
situation is similar to that encountered in frustrated 
magnetism^ where many different configurations have 
free energies that are, to a first approximation, identi- 
cal. In the present case, the degeneracy of the insulating 
states is lifted by taking fluctuations of the on-site parti- 
cle number into account. 

The validity of the approach we present is controlled 
by the size of the fluctuations in the superfluid order pa- 
rameter, which will be smaller in higher dimensions. At 
least in three spatial dimensions (3D), the fluctuations 
are expected to be relatively small, except close to the 
superfluid phase boundary. The approach presented here 
is therefore likely to be valid deep within both the insu- 
lating and superfluid phases, but necessarily breaks down 
in the vicinity of the phase transition to the superfluid. 
As usual, numerical studies are required to provide reli- 
able results for the exact positions of the various phase 
boundaries. 
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The model that we use includes only the lowest band in 
the optical lattice potential, which is appropriate for the 
temperatures and hopping strengths that are treated. An 
important conclusion of this work is that the maximum 
transition temperature to a magnetically ordered insula- 
tor is limited by the instability to superfluidity, and to 
occur for hopping strengths well within the regime where 
a one-band model is applicable. This contrasts with the 
fermionic case^ where higher bands must be taken into 
account. We furthermore restrict throughout to the spa- 
tially homogeneous case and ignore the effects of trap- 
ping. In the presence of an external parabolic trapping 
potential, the results presented here apply to the Mott 
insulating regions in the resulting tiered 'wedding cake' 
structured 

In Section [Til we introduce the model that is used 
throughout and review the limit where it can be de- 
scribed in terms of localized moments. In Section IIII1 
Landau theory is applied to characterize the possible 
phases and transitions that the model describes. In Sec- 
tion HVi we derive an approximation method based on 
a field-theoretical approach, and present the phase di- 
agrams that result. We conclude in Section [V] with a 
summary and some comments about experimental real- 
ization and detection of these phases. 



This model is equivalent to that studied by Soyler et al.^ 
and has the advantage of reducing the on-site Hilbcrt 
space to that of the corresponding problem for fermions. 

Our primary interest will be 'magnetic' phases, where 
the spin degrees of freedom order, and so we restrict to 
the case with mean filling of one particle per site and with 
no population imbalance. The reduction of the prob- 
lem to one of hard-core bosons introduces an additional 
particle-hole symmetry, and so the occupation number 
can be fixed with a chemical potential fi = U/2. 

It will sometimes be convenient to include a fictitious 
external 'magnetic field' h, which couples to the on-site 
spin through a Zceman term, 

H h = - hi ■ Si , (4) 

i 

and allows one to study the instability towards magnetic 
ordering. As we will show, this model can exhibit both 
ferromagnetic and (antifcrromagnetic) Neel order, and so 
hi will be taken as either uniform or staggered to describe 
these phases. 

B. Perturbation theory in t/U 



II. MODEL 

A. Hamiltonian 

We consider two species of bosons on a square or cubic 
lattice, described by the Hubbard model: 

H = - E *S4A« + E V^blblb n b la , (1) 

where sites are labeled by i and j, and species by a, 7 <E 
{1,2}. The hopping matrix element i™- is equal to t a > 
if i and j are nearest-neighbors and zero otherwise. 
Even for t\ ~ t^-, the Hamiltonian docs not in general 
have SU(2) symmetry, instead having only U(l) x U(l) 
symmetry under independent phase rotations for the two 
species of bosons, corresponding to conservation of both 
particle numbers separately. By analogy to the SU(2)- 
symmetric case, the 'spin' on site i can be defined in 
terms of the Pauli matrices <x as 

S i = biaVarrbvy ■ (2) 

The full spin-rotation symmetry is explicitly broken down 
to the set of rotations in the x-y plane. 

The calculations can be simplified considerably by re- 
stricting to the case where the intraspecies repulsion, 
V aa , is much larger than the interspecies repulsion, V ai 
for a/7. Taking the limit of infinite intraspecies repul- 
sion, we can describe both species by hard-core bosons 
and write 

n = - E wit* + u E *™ • ( 3 ) 

ij i 



In the limit t a /U <C 1, particle- number fluctuations 
are strongly suppressed and the physics is well described 
by an effective spin model. Restricting to the subspace 
where there is precisely one boson per site, the spin op- 
erator Si describes a moment of S = i, and an effective 
Hamiltonian can be derived using perturbation theory in 
t/U, leading to&i2 

Wspi„ = £ [-j x (s*s* + sfsf) + rsts*] . (5) 

(ij) 

The coupling constants J = At\ti/U and J z = 2{t\ + 
t^)/U are both positive,— and result from virtual tunnel- 
ing processes, where a particle hops back and forth be- 
tween two sites (J z ) or two particles with opposite spin 
on adjacent sites exchange places (J x ). 

The phase structure in this limit is therefore described 
by a quantum spin-^ XXZ model: At zero temperature, 
there is a first-order transition between a ferromagnet 
with (Si) in the x-y plane for J- 1 3> J z and a Neel state 
with (Si) along the z axis for J z ;§> J- 1 . Both give way 
to a paramagnetic phase at a critical temperature Tq 
proportional to the coupling. 

This perturbative calculation is only applicable in the 
limit of weak tunneling, and so cannot describe transi- 
tions to a supcrfluid state. It nonetheless provides the 
important insight that the tendency towards magnetic or- 
dering is caused by the enhancement of particle-number 
fluctuations between magnetically aligned sites. In Sec- 
tion |lVl we show how this physics can be incorporated in 
a calculation of the phase diagram that is valid beyond 
the perturbative limit. 
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III. LANDAU THEORY 

To understand the general form of the phase diagram, 
and the nature of the phase transitions that are possible 
between the different types of order, it is useful to con- 
sider the predictions of Landau theory. We first identify 
order parameters to distinguish the various phases of in- 
terest, then construct an expansion for the free energy in 
terms of these, including all terms allowed by symmetry. 

As noted above, the model is capable of both ferromag- 
netic order, where the spins are aligned in the x-y plane 
(and which will be referred to as the 'XY-ferromagnet'), 
and Ncel order, with the spins aligned with the z axis. We 
therefore define the uniform magnetization to ~ (Si) and 
the staggered magnetization n ~ ?/;{£f) (where r\i = ±1 
on the two sublattices) . Note that fixing the total parti- 
cle number of the two species to be equal (and assuming 
the absence of phase separation^) implies that m lies in 
the x-y plane, and we further assume that the Neel order 
is aligned with the z axis. 

As in the spinless case, the Hubbard model can also 
demonstrate supcrfluid phases, in which the phase ro- 
tation symmetry of one or both boson species is spon- 
taneously broken, and which are described by the two 
superfluid order parameters Xa ~ (bia)- The phase in 
which both are nonzero, which will be referred to simply 
as the superfluid, necessarily also has nonzero to, with 
its direction in space determined by the relative phase of 
Xi and X2- 

In their numerical study of the zero-temperature phase 
diagram of this model, Soyler et al.— also find phases 
with superfluidity of only one of the two species coex- 
isting with checkerboard density order. Such a 'checker- 
board superfluid' has xi ( sa y) an d n ^ 0, while X2 
and to vanish. 

We first treat the supcrfluid and ferromagnetic phases, 
for which the Landau theory involves the order parame- 
ters to, xi and X2j an d takes the general form 

£m, x = nxlxi + 7-2X2X2 + M n(XiXi) 2 + u 2 2(X2X2) 2 
+ ui2(XiXi)(X2X a ) + i?|TO| 2 + U(\m\ 2 ) 2 

+ D(m z ) 2 +"/?Ti- (x^x) + ■■ ■ ■ (6) 

The term with coefficient D > constrains to to lie in 
the x-y plane, the last term is the lowest-order coupling 
allowed between the magnetic and superfluid orders, and 
the ellipsis represents further terms of higher order. Ac- 
cording to Landau theory, the phase structure can be 
determined by minimizing £ m ,x a s a function of the pa- 
rameters appearing in this expression, which are in turn 
(unknown) functions of the physical parameters appear- 
ing in the Hamiltonian. 

In the superfluid phase, both xi and X2 become 
nonzero, with the relative magnitude determined by the 
ratios between the coefficients ri 2 and 1*11,22,12- Their 
relative phase remains arbitrary, but is fixed once the di- 
rection of m in spin-space is determined: for 7 > 0, to 
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FIG. 1: (Color online) Phase diagram predicted by Landau 
theory, determined by minimizing the free-energy function 
£m, x given in Eq. Q, as a function of A and a. (The pa- 
rameters B = b — —a = 1 are fixed without loss of gen- 
erality, and higher-order terms are set to zero.) The phases 
included are the spin-disordered Mott insulator (MI) , the XY- 
ferromagnetically ordered MI, labeled 'XY MP, and the su- 
perfluid, which also has XY-ferromagnetic order. The solid 
lines dividing the two MI phases and dividing the superfluid 
from the MI phases for A > j and A < — | are continuous 
transitions. The thick line shows a first-order transition be- 
tween the superfluid and MI phases, which is surrounded by 
a hysteretic region, shown shaded and surrounded by dashed 
lines. 



and x^&X are parallel, while for 7 < 0, they are antipar- 
allcl. The minimum of C„ 
the simplified expression 



allel. The minimum of C rn . x can therefore be found from 



a X 2 + bx + Am 2 + BrrT + amx 2 



(7) 



where the vector structure and complex phase of to and 
Xa have been eliminated by requiring that C m . x be min- 
imized: the minimum always occurs when the vectors to 
and x OX are aligned (and lie in the x-y plane). Note 
that a < 0, without loss of generality, since this mini- 
mization requirement fixes the relative orientation of to 
and x^X- 

The phase structure implied by this Landau action is 
shown in Figure U as a function of the parameters a 
and A, with the assumption that B and b are positive. 
The phase transitions separating the three phases are in- 
dicated with solid lines, which are thin for continuous 
transitions and thick for first-order. The shading indi- 
cates the region of hysteresis, where two local minima 
of C' m x exist and the stable phase is determined by the 
global minimum. 

In the plot, the coefficients b = B = —a = 1 have been 
fixed without loss of generality (assuming positive b and 
B; by rescaling to, x and the overall scale of £' mx ), and 
the coefficients of the remaining terms, represented by 
the ellipsis in Eq. ((?]), have been set to zero. As usual 
in Landau theory, the positions of the continuous transi- 
tions depend only on the coefficients of the first few terms 
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in the expansion, but the precise shape of the hysteretic 
region and the position of the first-order line also depend 
on the higher-order coefficients. (For B or b negative, 
further terms in the expansion become important and 
first-order boundaries will extend over a larger region of 
the diagram.) 

The following general conclusions of Landau theory are 
illustrated in Figure [T] Firstly, the three phases, spin- 
disordered Mott insulator (MI), XY-ferromagnetic MI, 
and supcrfluid, meet at a point. Secondly, the transition 
between the two insulating phases, across which m be- 
comes nonzero, is continuous (assuming that the quartic 
coefficient B is positive). Finally, the transitions from 
the insulating phases to the superfluid, across which xi 
and X2 become nonzero, are of first order in a region sur- 
rounding the point where all three phases meet, but can 
be continuous elsewhere. In Section llVl we will develop 
an approximate treatment of the microscopic model and 
show that its conclusions are in agreement with those 
presented here (see, in particular, Figure [3]). 

It is straightforward to include the Neel phase in the 
above analysis, by considering couplings of the order pa- 
rameter n to m and Xa- Due to the presence of the stag- 
gering factor rji, the only allowed terms involve n 2 , a 
scalar under the symmetry group. As a consequence, 
Landau theory predicts that a direct transition from a 
Neel phase, with n ^ and fh = 0, to either the su- 
pcrfluid or XY-ferromagnet must be of first order. This 
prediction is again in agreement with the results of Sec- 
tion [IV] (Figures H and [SJ . 

Finally, the checkerboard supcrfluid, with Xi 7^ and 
n / 0, can be connected by a continuous transition to 
the Neel phase, equivalent to the standard supcrfluid- 
insulator transition of a single species of boson. A direct 
transition to the uniform superfluid, where both xi an( i 
Xi are nonzero, is necessarily of first order. These obser- 
vations are consistent with the results of quantum Monte 
Carlo simulations^ 

IV. FIELD-THEORETICAL APPROACH 

As seen in Section IIIB1 perturbation theory in t a /U 
predicts a critical temperature Tq for magnetic ordering 
proportional to t 2 a jU . Maximizing Tq therefore requires 
increasing the hopping beyond the limit where this lead- 
ing order result is valid. While it is possible to continue 
the expansion in t a /U to higher order, generating fur- 
ther couplings between spins, such an approach cannot 
describe the superfluid phase. We instead use an approx- 
imation based on the mean-field theory that has been to 
applied the spinless case^ 

A. Partition function 

The thermodynamic properties are described in terms 
of the partition function, defined by Z = Tre _/3W , where 
the trace is over the full Hilbcrt space for all sites, and 
f3 — l/T is the inverse temperature. (Here and through- 



out, we set fee = 1.) The mean- field theory for the Bose- 
Hubbard model can be derived starting from a Hubbard- 
Stratonovich transformation) 7 ! 8 ' 14 which rewrites Z in 
terms of an integral over a complex field i\) a : 



J V 2 ip cxp — 


J 


f P 

' dT1pi(T)t- l 1p{T)+Sl[i)] 



J V 2 ijj exp- 




J 








(8) 



The auxiliary field ip a has the same vector structure as 
the boson operator bi a , and we use the shorthand 

^( r )t-V(r) =S4W(t; 1 W; a (r). (9) 

ij,ct 

Once the hopping term has been decoupled using 
the Hubbard-Stratonovich transformation, the remaining 
terms in the Hamiltonian each act only at a single site. 
The effective action 51 can therefore be found in terms of 
the solution of a one-site problem: 

f2[V>] = -5^togTkT T e-tf dT ««W«<»«) , (10) 

i 

where T T denotes ordering in imaginary time and 

Ui{ip a ) = Unan i2 - n(n a + n i2 ) - K ■ &j Q o ! a , 7 6 i7 

-J2^ha+Ma) (11) 

is the local effective (time-dependent) Hamiltonian. 

1. Gaussian approximation 

To reduce the functional integral in the numerator of 
Eq. ([5]) to a tractable form, we approximate the effective 
action Q[ip] by an expansion up to quadratic order around 
its minimum. The validity of this approximation is con- 
trolled by the size of the fluctuations, which are expected 
to be larger closer to the superfluid transition and in 
fewer spatial dimensions. In particular, in the supcrfluid 
phase in 2D at nonzero temperature, the superfluid order 
parameter is completely eliminated by fluctuations, and 
this approach is not expected to be applicable. In 3D, 
we expect the approximation to give reasonable quanti- 
tative results, except in the region close to the transition 
to the superfluid, where the gap to single-particle excita- 
tions vanishes and a quadratic approximation ceases to 
be valid. This approach also omits fluctuations of the 
magnetization, which become large in the region close to 
the magnetic ordering transition. 

For the purposes of studying the magnetic ordering in 
the insulator and finding the phase boundary to the su- 
pcrfluid, it is sufficient to expand f2 around the point 
where ip a = 0. Using time-dependent perturbation the- 
ory, Eq. pop can be expanded in powers of ip as 
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n = -^logTrc-^W + ^EE^ G ^M^ + W 4 , (12) 



uj 2,07 

whcr o 16 i 17 

X- e -^„ \^ { (n\bl\n')(n'\b^\n) ^ (n|6 7 K 



ilu + e n > - e n — vjj + e n > - e„ 



G iai {iuj) = — - ^ n _ ptn (13) 



and denotes a sum over all (bosonic) Matusbara 
frequencies uj. In Eq. (| 13[) . the indices n and n' la- 
bel eigenstates of the single-site Hamiltonian Wj(0), and 
e n and e n ' are the corresponding eigenvalues. In the 
zero-temperature limit, the sum over n' becomes a sum 
over excitations, in this case double and zero occupation, 
above the on-site ground state. (The restriction of the 
on-site Hilbert space to that of hard-core bosons signifi- 
cantly simplifies the calculation of this quantity.) 



Once the expansion in Eq. (|12p has been truncated 
to quadratic order, both integrals over ip in Eq. ([5]) are 
Gaussian and can be calculated to give the free energy, 
T = - A log Z. We find T = JF + T x , with 



lim| z |_, c 



= 1, 



Q(iw) 



1 

2/3 



^log(cosh/?2;p — 1) 



^log(cosh/3z Q - 1) 



(16) 



where zp and zq are the zeros of P(z) and Q(z). 

The calculation of the free energy T , as a function of 
the applied field, the temperature and the parameters in 
the Hamiltonian, is therefore reduced to the problem of 
calculating the matrix G and integrating over k. 



^H-^logTre-^ ), ( 14 ) 

^iW = ^E lo g det [ 1+tG ( iw )] . ( 15 ) 



where t and G are matrices in sites and flavor indices. 
Note that the dependence on the external field h is due to 
the Zeeman term in Jii and the resulting dependence on h 
of the eigenstates and eigenvalues appearing in Eq. (|13[) . 

The saddle-point contribution to the free energy, J-q, 
has no dependence on the hopping t a and only contains 
information about the on-site state. The Gaussian fluc- 
tuations about this saddle point give T\ , which takes into 
account the particle-number fluctuations that are crucial 
for describing the magnetic ordering transitions. In the 
case where hi is uniform, T\ can be written as a sum 
over eigenvectors of t, labeled by momentum k; in the 
thermodynamic limit this becomes a c?-dimensional inte- 
gral over k. When the applied field is staggered (such as 
when considering Nccl ordering), there is mixing between 
momenta separated by the appropriate reciprocal lattice 
vector, but the same approach can be used. 

The Matsubara sum in Eq. (fT5"|) converges because the 
determinant tends (sufficiently rapidly) to 1 for large fre- 
quencies. Since G(z) and hence det[l + tG(z)] arc ra- 
tional functions, the sum can be calculated using the 
result that for any polynomials P and Q such that 



2. Single-particle gap 

The perturbative calculation in Section Ul Bl shows that 
magnetic ordering results from the enhancement of hop- 
ping between sites with aligned spins. One therefore ex- 
pects that the gap A to single-particle excitations should 
depend on the magnetization, and that capturing this 
effect is crucial to describing the magnetic phases. By 
inverting the Hubbard-Stratonovich transformation, the 
propagator of the bosons 6; Q can be related^ to that 
of the field ip a , and is given within our approximation 
scheme by [G(icj) _1 + t] _1 . The single-particle excita- 
tions are described by the poles of this expression for 
real z = iu>, and the gap A is given by the smallest exci- 
tation energy. Figure [5] shows a plot of this quantity in 
the insulating phase, in the presence of a uniform applied 
field in the x direction. 

For the temperatures shown in the plot, thermal 
particle-number fluctuations arc strongly suppressed, 
and the lowest-energy on-site configuration has a sin- 
gle particle. For h x <C T, this particle can have any 
spin orientation, while for h x 3> T, the spin is fixed to 
point along the x axis. For h x considerably smaller than 
T, increasing the field makes the on-site configuration 
more spin-polarized, allowing for more coherent particle- 
number fluctuations, and decreasing the gap. This con- 
tinues until the applied field is somewhat larger than 
the temperature and the on-site state is essentially com- 
pletely spin polarized. The gap then increases with in- 
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FIG. 2: Single-particle gap A, in units of the on-site inter- 
action strength U, as a function of (uniform) applied mag- 
netic field h x , for three different temperatures. The hopping 
strengths are ti = ti = 0.04(7 and the lattice coordination 
is z — 6. For h x < T, the gap decreases rapidly with h x as 
the spins become ordered, enhancing particle-number fluctu- 
ations between neighboring sites. This fluctuation effect can 
lead to spontaneous magnetic ordering in the absence of a 
field and is captured by including Gaussian corrections in the 
calculations of Section IIVI For h x 3> T, the magnetization 
is saturated and the Zeeman term causes the gap to increase 
linearly with \h x \. 



creasing field, since the on-site energy cost of particle- 
number fluctuations grows linearly with \h\. 



For vanishing applied field, T\ becomes 




for a uniform field applied at an angle 9 from the z axis. 
The free energy is minimized by 9 = ?, i.e., by taking the 
applied field to lie in the x-y plane, in agreement with the 
results of perturbation theory in Section Hi Bl For 6 = 
the free energy reduces to 

(18) 

The corresponding expression for a staggered field is min- 
imized by a field along the z axis, 9 = 0, for which the 
free energy is given by 

+y^-Atij . (19) 

These expressions for T agree with those given by 
Altman et al.,— in particular their Eqs. (33) and (41). 
Expanding in powers of the hopping^ one reproduces 
the results of perturbation theory, given above in Section 



2. Nonzero temperature 



B. Magnetic ordering 

The phase boundaries for magnetic ordering can be 
found using the expression for the free energy F(h) given 
in Eqs. 1141 and 1151 Different approaches are required for 
zero and nonzero temperatures. 



1, Zero temperature 



For j3 < co, G is a continuous function of the applied 
field, and the previous approach no longer applies. To de- 
termine where magnetic ordering takes place, we find the 
effective potential $ as a function of the magnetization 

m. First, we define h{rfi) as the solution of 



dT_ 

dh 



(20) 



h(rh) 



The effective potential $(m) is then given by the Legen- 
dre transform of J-(h): 



In the limit /3 — * co, the only term that contributes 
to the sum over n in Eq. (fT3"|) is where \n) is the ground 
state of the on-site Hamiltonian Hi(0). There is therefore 
a discontinuity at zero applied field, because 7Yi (0) then 
has a two-fold degenerate ground state. This degeneracy 
results from the spin degree of freedom, and is split by an 
infinitesimal h. The limit of the free energy T as h —> 
therefore depends on the direction in spin-space, and on 
the configuration of hi in real space. The ground state is 
found by minimizing T with respect to the direction of 
h in the limit that its magnitude goes to zero. 



$(m) = m • h(m) + J"(/i(m)) . 



(21) 



The expression for the free energy, J- = J- + J-i, can 
be viewed as the first two terms in an expansion in the 
order of the fluctuations, and we similarly seek $(m) as a 
series in fluctuations. For this we require only the zero- 
order expression ho, satisfying —^(hoijn)) = m: the 
expression for $(m) is given by 

$o(m) + $i(m) = m- h (m) + T {h (m)) + ^i(K (m)) . 

(22) 
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This result for <f> differs from one that would be ob- 
tained by performing the Legendre transformation di- 
rectly with T = J-q + J- 1 , and is equivalent to an RPA- 
like summation. For example, (the exact result for) the 
quadratic coefficient in a series expansion of $ is given by 
$"(0) = -[^"(O)] -1 = x" 1 , the reciprocal of the mag- 
netic susceptibility. The expression in Eq. (|22p . which 
gives 4>"(m) consistently in powers of the fluctuations, is 
therefore equivalent to rcsumming a geometric progres- 
sion: 



$"(0) = (Xo + Xi)~ 



1_ xi 



Xo 



Xi 
Xo 



(23) 



As in the familiar case of the Stoner criterion, the phase 
transition to a magnetic state occurs when the suscepti- 
bility diverges, and the 'correction' has the same magni- 
tude as the leading order result. 

To find the magnetic phase structure within the Mott 
insulator, we therefore calculate <& according to Eq. (|22|) 
and minimize it with respect to the magnetization to. 
If the minimum occurs for zero magnetization, then the 
system is magnetically disordered. If it occurs for uni- 
form nonzero to, the system is ferromagnetically ordered, 
while if it occurs for staggered to, the system is antifer- 
romagnctically ordered. 



C. Superfluid— insulator transitions 
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FIG. 3: Phase diagram on the simple cubic lattice as a func- 
tion of temperature T and hopping t — t\ — t-z- The phases 
are labeled as in Figure [1] with 'XY MI' referring to a fer- 
romagnetically ordered Mott insulator (MI) where the mag- 
netization vector lies in the x-y plane. The two insulating 
phases are separated by a continuous phase transition (thin 
line), while the superfluid-insulator transitions are either con- 
tinuous or first order (thick line). The structure of the phase 
diagram, including the order of the transitions, agrees with 
the predictions of Landau theory shown in Figure [1] (At 
least in this approximation, the transition between the spin- 
disordered MI and superfluid remains first order for arbitrarily 
large t, indicating that the point A = j in the Landau theory 
is never reached.) 



As described in Section IIV A 21 the energies of single- 
particle excitations can be determined by locating the 
poles of the bosonic propagator. At the transition from 
the insulator to the superfluid, the gap vanishes and the 
Gaussian approximation scheme breaks down. As in the 
spinlcss case,— it is nonetheless possible to estimate 
the position of the phase boundary as the point where 
the lowest-order estimate for the gap vanishes. 

In the present case, however, it should be noted that 
the single-particle gap is dependent on the magnetiza- 
tion, as illustrated in Figure [H and it is therefore possi- 
ble for the gap to vanish first at a value of magnetization 
away from the minimum of the effective action. In the 
region where the two hopping strengths t\ and t2 are 
of similar magnitude, and ferromagnetic ordering is fa- 
vored, the gap is in fact found to decrease with increasing 
magnetization, as expected from the analysis of Section 
Mil As predicted by Landau theory, this can lead to a 
first-order transition, where the appearance of a super- 
fluid order parameter is accompanied by a discontinuous 
increase in the magnetization. 

D. Phase diagrams 

This approximate treatment allows the phase structure 
to be found, as a function of the two hopping strengths, 
t\ and tii and the temperature T, which we measure in 



units of the interspecies repulsion U . Sections through 
the three-dimensional parameter space are shown in Fig- 
ures [31111 and[5l where, as in Figure [U first-order transi- 
tions are indicated by thick lines, and continuous transi- 
tions by thin lines. All three phase diagrams have been 
calculated assuming nearest-neighbor hopping on a sim- 
ple cubic lattice. 

Figure[3]shows the phase diagram in the plane of equal 
hopping for the two species, t\ = t2- As in the pcrtur- 
bative calculation of Section [LTbI ferromagnetic order is 
preferred over Nccl, and the phase diagram contains the 
XY-ferromagnet, as well as the spin-disordered Mott in- 
sulator and the superfluid. The general phase structure is 
therefore directly comparable to, and in agreement with, 
that predicted by Landau theory shown in Figure [TJ 

It should be noted that the maximum temperature for 
the XY-ferromagnet is limited by the presence of super- 
fluid order. While the transition temperature for mag- 
netic order increases with hopping, as in the perturbative 
calculation of Section Ul Bl that for superfluidity increases 
more rapidly. The superfluid replaces the ferromagnetic 
insulator entirely for moderate t/U, at which point the 
global maximum of the ferromagnetic transition temper- 
ature is reached. 

Figure [4] shows a phase diagram with fixed t\ + 1% = 
0.08U, as a function of hopping t\ and temperature T. 
Near t\ = 0.04L 7 , where the hopping strengths are equal, 
the XY-ferromagnet has lower free energy, but for larger 
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FIG. 4: Phase diagram as a function of temperature T and 
hopping ti, with fixed ti + ti = 0.08(7. For this value of the 
total hopping, there is no superfiuid phase. For high temper- 
atures, the system is a Mott insulator (MI) without spin or- 
der. For lower temperatures and ti ~ ti, a ferromagnetically 
ordered state (labeled XY MI) is favored, with the magneti- 
zation vector in the x-y plane. For more imbalanced hopping 
values, a Neel state is preferred, with a two-sublattice order 
and the staggered magnetization vector parallel to the z axis. 
As expected from the Landau theory of Section Mil the tran- 
sition from the spin-disordered insulator to either magnetic 
phase is continuous, but the transition between the two mag- 
netic orderings is of first order. 



disparity in the hoppings, the Neel phase is instead fa- 
vored. As predicted in Section [Till there is a first-order 
transition between the two, at which both magnetic or- 
der parameters change discontinuously. The value of t\ 
at this transition line is only weakly dependent on tem- 
perature, with the extent of the Neel phase increasing 
slightly with T. The two magnetically ordered phases 
survive only until a temperature of order 1% of U, where 
there is a continuous transition into the magnetically dis- 
ordered phase. (An attempt to raise the critical temper- 
ature further by increasing the hopping leads instead to 
the superfiuid phase.) 

In Figure the phase diagram is shown for fixed tem- 
perature T = 0.01Z7, with all four phases appearing for 
suitable choices of the hopping strengths t\ and ti. At 
this temperature, the transition from the ferromagnetic 
insulator to the superfiuid is continuous, but, as pre- 
dicted in Section [iTTl the transition from the Neel phase 
to the superfiuid is of first order. 



V. DISCUSSION 

This work has presented a study of a system of two 
species of bosons moving in a lattice, based on two 
approaches. First, the standard techniques of Landau 
theory were used to determine the nature of the possi- 
ble phase boundaries. Then an analytic framework was 
introduced, based on a field-theoretical representation 
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FIG. 5: Phase diagram at fixed temperature T = 0.01/7, as a 
function of the two hoppings ti and t2- At this temperature, 
all four phases appear for certain values of the hoppings, and 
the transition between the XY-ferromagnet and the super- 
fluid is continuous. As in Figure [4] an antiferromagnetically 
ordered Neel state is favored for more imbalanced hoppings. 
The transition from this phase to the superfiuid is of first or- 
der, in agreement with the conclusions of Landau theory. This 
phase diagram is directly comparable to Figure 5 of Ref. HcL 
which is plotted at zero temperature, where there is no spin- 
disordered MI. 



of the partition function. This builds upon the stan- 
dard mean-field theory for the spinless Bose-Hubbard 
model,— & but takes into account the important effects 
of fluctuations, leading to magnetic order. Extensions to 
the method could be used to describe other bosonic mod- 
els and to calculate other properties, such as dynamic 
correlations. 

An important conclusion of this work is that the maxi- 
mum temperature at which magnetically ordered insulat- 
ing phases can survive is determined by the presence of 
the superfiuid phase. While the critical temperature for 
magnetic ordering increases with hopping, in line with 
the predictions of a simple perturbative calculation, the 
same is true of the superfiuid transition temperature. For 
moderate values of the hopping, reducing the tempera- 
ture from the spin-disordered Mott insulator leads di- 
rectly to the superfiuid phase, without passing through 
a magnetically ordered insulator (sec Figure [3]). 

This implies that an experimental realization of the 
magnetically ordered insulators will require temperatures 
in the lattice on the order of 1% of the on-site repulsion 
U. At these temperatures, thermal particle-number fluc- 
tuations arc completely frozen out, and the only entropy 
results from spin fluctuations. Reaching these phases will 
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therefore likely require a considerable reduction in the 
entropy after loading into the lattice, using techniques 
such as those proposed by Bcrnier et al.— To give precise 
predictions for the temperature and entropy constraints, 
further numerical work is needed, and the effects of the 
trapping potential, which has been neglected here, should 
be taken into account. 

As pointed out by Altman et al. ) 10 i 18 the magnetically 
ordered phases could be detected by analyzing correla- 
tions of the shot noise in time-of- flight measurements. It 
would also be interesting to use recently developed tech- 



niques of photocmission spectroscopy for cold atomsi^ to 
investigate the effect of the magnetic ordering transition 
on the single-particle Green function.— 
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The approximation scheme used in this work is not capa- 
ble of describing such phases and further numerical work 
is required to determine their stability to thermal fluctua- 
tions. 

This should be contrasted with the fermionic case, where 
an extra minus sign from the exchange of fermions leads 
to antiferromagnetic J x as well as J z . 
Infinite intraspecies repulsion means that real-space segre- 
gation of the two species is likely to be extremely unfavor- 
able energetically. 



